# -*- coding: utf-8 -*-
# Name        модуль1
# Purpose
#
# Author      daimoon
#
# Created     02.11.2013
# CopYright   (c) daimoon 2013
# Licence     <Your licence>
#-------------------------------------------------------------------------------
import sys
sys.path.append("./modules")
import numpy as np
from math import *
import const

class upper_estimate(object):
    def __init__(self, n_ob, n_el, l, l1, h, b, phi, c, gam, q, k_l, k_l1):
        self.n_ob = n_ob
        self.n_el = n_el
        self.b = b
        self.l = l
        self.l1 = l1
        self.k_l = k_l
        self.k_l1 = k_l1
        self.h = h
        self.gam = gam
        self.q = q

        self.u = np.zeros(n_el)
        self.v = np.zeros(n_el)
        self.M = np.zeros(n_el)
        self.Y = np.zeros(n_el)
        self.X = np.zeros(n_el)
        self.uv = np.zeros(n_el)
        self.Mq = 0
        self.Msw = 0
        self.i = 1
        self.j = 1

        self.A = (1 + sin(phi)) / (1 - sin(phi))
        self.Cb = 2 * c * cos(phi) / (1 - sin(phi))
        self.s2 = self.sqr((self.A + 1) / (self.A - 1))
        self.alf = pi / 2 - atan(sqrt(self.s2 - 1))

    def set_2(self, l, n):
        if n < n * 0.5:
            n = -n
        l = l + n * l
        return l

#    def set(self, l, n):
#        if random.random()>0.4:
#            l_new = self.set_2(l, n)
#        else:
#            l_new = self.set_1(l, n)
#        if l > 0.1 * l_new:
#            l = l_new
#        return l

    def set(self, l, n):
        l_new = self.set_2(l, n)
        if l > 0.1 * l_new:
            l = l_new
        return l


    def sqr(self, x):
        sqrx = x * x
        return sqrx

    def block_1(self):
        self.u[1] = 1
        self.v[1] = 0.5 * self.A * self.b * self.u[1] / self.h
        self.M[1] = self.Cb * self.b * self.u[1] / 4
        self.Y[1] = 0
        self.X[1] = self.h
        self.Y[2] = 0.5 * self.b
        self.X[2] = 0
        self.i = 1
        self.j = 2

    def met1(self):
        if self.j > 2:
            self.j += 1
        self.uv[self.i] = sqrt(self.u[self.i] * self.u[self.i] + self.v[self.i] * self.v[self.i])
        de = atan(self.u[self.i] / self.v[self.i])
        self.Y[self.j+1] = self.Y[self.j-1] + self.l1 * cos(self.alf + de)
        self.X[self.j+1] = self.X[self.j-1] + self.l1 * sin(self.alf + de)
        self.M[self.j] = 0.5 * self.Cb * abs(self.uv[self.i]) * self.l1 * (1 - cos(pi * 0.5 - self.alf))

    def met3(self):
        ll = sqrt(self.sqr(self.X[self.j+1] - self.X[self.j]) + self.sqr(self.Y[self.j+1] - self.Y[self.j]))
    ##   блок 3
    ##   скорости в доп-й системе
        ts = self.X[self.j+1] / ll
        tc = (self.Y[self.j] - self.Y[self.j+1]) / ll
        vd = self.v[self.i] * tc - self.u[self.i] * ts
        ud = self.u[self.i] * tc + self.v[self.i] * ts
        eY = vd / ll
        gXY = 2 * ud / ll
        eX = (sqrt(self.sqr(self.s2 + 1) * self.sqr(eY) - (self.s2 - 1) * ((self.s2 - 1) * self.sqr(eY) - self.sqr(gXY))) - (self.s2 + 1) * eY) / (self.s2 - 1)
        ebs3 = (eX + eY) / 2 - 0.5 * sqrt(self.sqr(eX - eY) + self.sqr(gXY))
        self.X[self.j+2] = 0
        self.Y[self.j+2] = self.Y[self.j] + self.l
        XX = self.l * ts
        YY = ll + self.l * tc
    ##  скорости в следующей верши­е тр-к  в доп-й системе
        self.i = self.i + 1
        self.j = self.j + 1
        ud = eX * XX + gXY * YY / 2
        vd = gXY * XX / 2 + eY * YY
    ##  скорости в следующей верши­е тр-к  в глоб-й системе
        self.u[self.i] = ud * tc - vd * ts
        self.v[self.i] = vd * tc + ud * ts
    ##  определе­ие площ ди тр-к
        ll1 = sqrt(self.sqr(self.X[self.j+1] - self.X[self.j]) + self.sqr(self.Y[self.j+1] - self.Y[self.j]))
        p = 0.5 * (ll1 + self.l + ll)
        S = (p * (p - ll1) * (p - self.l) * (p - ll))
    ##  мощ­ость в ­ечет­ом блоке
        self.M[self.j] = 2 * self.Cb * (-ebs3) * S
        #################################################################################*2

    def met4(self):
        de = atan(self.u[self.i] / self.v[self.i])
        self.Y[self.j+1] = self.Y[self.j-1] + self.l1 * cos(self.alf + de)
        self.X[self.j+1] = self.X[self.j-1] + self.l1 * sin(self.alf + de)
        self.Y[self.j+1] = self.Y[self.j-1] - self.X[self.j-1] * (self.Y[self.j+1] - self.Y[self.j-1]) / (self.X[self.j+1] - self.X[self.j-1])
        self.X[self.j+1] = 0
        self.l1 = sqrt(self.X[self.j-1] * self.X[self.j-1] + self.sqr(self.Y[self.j+1] - self.Y[self.j-1]))
        k1 = self.j + 1
        for k in range(1, self.j, 1):

            if self.Y[self.j+1] < self.Y[k+1]:
                self.j = k + 1
                jtr = self.j - 2
                ij = round(self.j / 2)
                k2 = round(self.j / 2) - 1
                self.X[self.j] = self.X[k1]
                self.Y[self.j] = self.Y[k1]
                self.l1 = sqrt(self.sqr(self.X[self.j-1]) + self.sqr(self.Y[self.j] - self.Y[self.j-1]))
                self.M[self.j-1] = 0.5 * self.Cb * abs(self.uv[k2]) * cos(pi / 2 - self.alf) * self.l1
                self.u[ij] = 0
                lpos = self.Y[self.j] - self.Y[self.j-2]
                sumu = (self.u[1] + self.u[ij-1]) / 2

                for ih in range(2, ij-2, 1):
                    sumu = sumu + self.u[ih]

                self.Mq = self.q * self.l * sumu + self.q * self.u[ij-1] * lpos / 2
                jh=0
                usr = np.zeros(ij * 2)
                st = np.zeros(ij * 2)

                for ih in range(1, ij-1, 1):
                    jh = jh + 1
                    usr[jh] = self.u[ih] / 3
                    jh = jh + 1
                    usr[jh] = (self.u[ih] + self.u[ih+1]) / 3
                    usr[jh] = self.u[ij] / 3
                    self.Msw = 0

                for ih in range(1, jtr):
                    o1 = sqrt(self.sqr(self.X[ih] - self.X[ih+1]) + self.sqr(self.Y[ih] - self.Y[ih+1]))
                    o2 = sqrt(self.sqr(self.X[ih+2] - self.X[ih+1]) + self.sqr(self.Y[ih+2] - self.Y[ih+1]))
                    o3 = sqrt(self.sqr(self.X[ih] - self.X[ih+2]) + self.sqr(self.Y[ih] - self.Y[ih+2]))
                    p = 0.5 * (o1 + o2 + o3)
                    st[ih] = sqrt(abs((p * (p - o1) * (p - o2) * (p - o3))))
                    self.Msw = self.Msw + st[ih] * self.gam * usr[ih] / 1000

                self.M[self.j] = 0.5 * self.Cb * abs(self.uv[self.i]) * cos(pi / 2 - self.alf) * self.l1
                self.j = self.j + 1
                jtr = self.j - 2
                ij =round((self.j - 1) / 2)
                self.u[ij+1] = 0
                lpos = self.Y[self.j] - self.Y[self.j-1]
                sumu = (self.u[1] + self.u[ij]) / 2

                for ih in range(2, ij-1, 1):
                    sumu = sumu + self.u[ih]

                self.Mq = self.q * self.l * sumu + self.q * self.u[ij] * lpos / 2
                jh = 1

                for ih in range(1, ij-1, 1):
                    usr[jh] = self.u[ih] / 3
                    jh = jh + 1
                    usr[jh] = (self.u[ih] + self.u[ih+1]) / 3
                    jh = jh + 1

                usr[jh] = self.u[ij] / 3
                self.Msw = 0

                for ih in range(1, jtr):
                    o1 = sqrt(self.sqr(self.X[ih] - self.X[ih+1]) + self.sqr(self.Y[ih] - self.Y[ih+1]))
                    o2 = sqrt(self.sqr(self.X[ih+2] - self.X[ih+1]) + self.sqr(self.Y[ih+2] - self.Y[ih+1]))
                    o3 = sqrt(self.sqr(self.X[ih] - self.X[ih+2]) + self.sqr(self.Y[ih] - self.Y[ih+2]))
                    p = 0.5 * (o1 + o2 + o3)
                    st[ih] = sqrt(abs((p * (p - o1) * (p - o2) * (p - o3))))
                    self.Msw = self.Msw + st[ih] * self.gam * usr[ih] / 1000

    def met5(self):
#        ugu = (self.u[1] - self.u[2]) / self.l
#        Mar=2 * sigt * af * hf * ugu
        Mm = 2 * (-self.Mq - self.Msw)

        for k in range(1, self.j - 1):
            Mm = Mm + self.M[k]
        Nwpr = Mm / self.u[1]
        Npriw = Nwpr * 100 / (self.b * self.b * self.gam)
        return Npriw

    def met6(self):
        u1 = np.zeros(self.n_el)
        v1 = np.zeros(self.n_el)

        for ih in range(1, self.j):
            self.i = round(ih / 2 + 0.1)
            k = round(ih / 2 - 0.1)
            if (self.i - k) == 0:
                v1[ih] = self.v[self.i]
                u1[ih] = self.u[self.i]
            else:
                v1[ih] = 0
                u1[ih] = 0
        return (v1, u1)

    def proc(self):
        """
        вход - self.l, self.l1, self.h, self.b, phi, c, self.gam, self.q
        выход - (Npriw, [self.X, self.Y, self.u, self.v, u1, v1])
        """
        self.block_1()
        p = 1
        n = 0
        end = 0
        while end == 0:
            if p == 1:
                self.met1()
                if (self.X[self.j+1] > 0):
                    self.met3()
                    self.l = self.set(self.l, self.k_l)
                    self.l1 = self.set(self.l1, self.k_l1)
                else:
                    p = 0
                    end = 1
            n += 1
            if n > self.n_ob: end = 1
        if p == 0 and n > 3:
            self.met4()
            Npriw = self.met5()
            v1, u1 = self.met6()
            X, Y, u, v = self.X, self.Y, self.u, self.v
        else:
            Npriw, X, Y, u, v, u1, v1 = const.bed_fitness, 0, 0, 0, 0, 0, 0
        return Npriw, X, Y, u, v, u1, v1

def estimate(l, l1, h, b, phi, c, gam, q, n_el1, n_el2, k_l, k_l1):
    est = upper_estimate(n_el1, n_el2, l, l1, h, b, phi, c, gam, q, k_l, k_l1)
    Npriw, X, Y, u, v, u1, v1 = est.proc()
    return Npriw, X, Y, u, v, u1, v1

def chi(out):
    m_n_x = np.eye(out[0].shape[0], 2)
    m_n_y = np.eye(out[0].shape[0], 2)
    j = -1
    for i in range(0, m_n_x.shape[0]):
        if out[0][i] != 0:
            j = j + 1
            m_n_x[j] = out[0][i]
            m_n_y[j] = out[1][i]
    m_n_x[j + 1] = min(out[0])
    m_n_y[j + 1] = max(out[1])
    return m_n_x, m_n_y


if __name__ == '__main__':
    h, l, l1 = 0.43, 0.52, 1.18
    b, phi, c, gam, q = 1.5, 30, 21, 18, 27
    phi = phi * (pi / 180)
    Npriw, X, Y, u, v, u1, v1 = estimate(l, l1, h, b, phi, c, gam, q, 30, 100, 0, 0)
    print(Npriw)